%
% Defining isotropic compliance and stiffness tensors using 
% complianceTensor  and  stiffnessTensor introduced in MTEX 5
%
% David Mainprice 9/03/2018
%
% Data from
% J.Fortin, Y.Gueguen, and A.Schubnel (2007)
% Effects of pore collapse and grain crushing on ultrasonic velocities and
% Vp/Vs, J. Geophys. Res., 112, B08207, doi:10.1029/2005JB004005.
%

% sandstone constants
Eo = 42.13; % Youngs Modulus in GPa
Ko = 21.3;  % Bulk Modulus in GPa
Go = 18.0;  % Shear Modulus in GPa
vo=(3*Ko-2*Go)/(2*(3*Ko+Go)); % Poissons ratio 

%% Isotropic tensor compliance Sij defined by Eo and vo

S11=(1/Eo); S12=(-vo/Eo); S44=2*(S11-S12);
% define matrix M
 M =... 
 [[  S11     S12    S12    0.0     0.0    0.0];...
 [   S12     S11    S12    0.0     0.0    0.0];...
 [   S12     S12    S11    0.0     0.0    0.0];...
 [   0.0     0.0    0.0    S44     0.0    0.0];...
 [   0.0     0.0    0.0    0.0     S44    0.0];...
 [   0.0     0.0    0.0    0.0     0.0    S44]];

CS_Tensor_Iso = crystalSymmetry('1',[1 1 1],...
 [ 90.0 90.0 90.0]*degree,'X||a*','Y||b',...
 'mineral','Isotropic tensor')

S_iso = complianceTensor(M,CS_Tensor_Iso)

%% Isotropic tensor stiffness Cij defined by Ko and Go
C11 = Ko+(4/3)*Go ; C12=C11-2*Go; C44=(C11-C12)/2;
% define matrix M
 M =... 
 [[  C11     C12    C12    0.0     0.0    0.0];...
 [   C12     C11    C12    0.0     0.0    0.0];...
 [   C12     C12    C11    0.0     0.0    0.0];...
 [   0.0     0.0    0.0    C44     0.0    0.0];...
 [   0.0     0.0    0.0    0.0     C44    0.0];...
 [   0.0     0.0    0.0    0.0     0.0    C44]];
CS_Tensor_Iso = crystalSymmetry('1',[1 1 1],...
 [ 90.0 90.0 90.0]*degree,'X||a*','Y||b',...
 'mineral','Isotropic tensor')
C_iso_solid = stiffnessTensor(M,CS_Tensor_Iso)
% Calculate  Compliance tensor by inversion
S_iso_solid = inv(C_iso_solid)
